

##################################################################################################################################################################################
#
#  Martini & Walter: "Learning from Precedent" - Plots      
#
##################################################################################################################################################################################


# Install/load required packages and set options:

require(ggplot2)
require(zoo)

options("scipen"=999, digits=4)

# Note: For R and package versions used, see end of script:


# Auto-detect and set working directory: 

DIR <- strsplit(rstudioapi::getActiveDocumentContext()$path, "/")[[1]]
DIR <- paste(DIR[1:(length(DIR)-1)], collapse = "/")
setwd(DIR)


##################################################################################################################################################################################


### PLOTS I: RHS measures - XR etc. => How well is Brexit going?


dev.off() # to reset plotting device, if needed
default_params <- par()
par(mgp=c(2.7,1,0)) #default: par(mgp=c(3,1,0))
label_size <- .7


Print <- F


{

# 1. Pound-Euro Exchange Rate

setwd(paste0(DIR,"/Data"))
files <- list.files(getwd(), ".txt")
XR_EUR <- read.table(files[grep("EUR_2022", files)], header = T, stringsAsFactors = F)
colnames(XR_EUR)[2] <- "eur_XR"

XR_EUR <- subset(XR_EUR, Date >= as.Date("2014-01-01")) 
XR_EUR <- subset(XR_EUR, Date < as.Date("2022-01-01")) 

setwd(paste0(DIR,"/Plots"))


if (Print == T) jpeg("rhs_1_XR.jpeg", units = "cm", width = 17, height = 16, res = 400)

plot(as.Date(XR_EUR$Date), XR_EUR$eur_XR, type = "l", col = "darkblue", ylim = c(1.07, 1.47), 
     xlab = "Timeline", ylab = "Spot exhange rate into Pound Sterling", main = "Pound Exchange Rate against the Euro") # blue

abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win elections => referendum coming
abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit

pos <- 1.47 
text(y = pos, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
text(y = pos, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 

if (Print == T) dev.off()

}


##################################################################################################################################################################################


# 2. UKIP and Tories polling vote shares

{
  
setwd(paste0(DIR,"/Data"))
library(lubridate)

ukip <- read.table("UKIP_polling_vote_share_monthly.txt", header = T, stringsAsFactors = F)
colnames(ukip)[1] <- "Month"
ukip <- ukip[-grep("2022", ukip$Month),] 
tories <- read.table("Tories_polling_vote_share_monthly.txt", header = T, stringsAsFactors = F)
colnames(tories)[2] <- "Tories_vote"
colnames(tories)[1] <- "Month"
tories <- tories[-grep("2022", tories$Month),] 

setwd(paste0(DIR,"/Plots"))


if (Print == T) jpeg("rhs_2_Ukip_Tories.jpeg", units = "cm", width = 17, height = 16, res = 400)

plot(as.Date(paste0(tories$Month, "-01")), tories$Tories_vote, type = "l", col = "darkblue", ylim = c(.5,56),
     xlab = "Timeline", ylab = "Polling results (%)", main = "UKIP and Conseratives Popularity", lwd = 1.75) # blue
lines(as.Date(paste0(ukip$Month, "-01")), ukip$Ukip_vote, type = "l", col = "steelblue2", lwd = 1.75) 
lines(as.Date(paste0(ukip$Month, "-01")), rowSums(cbind(tories$Tories_vote, ukip$Ukip_vote)), type = "l", col = "purple", lty = 2, lwd = 1.75) 

abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win elections => referendum coming
abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit

pos <- 56
text(y = pos, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
text(y = pos, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 

if (Print == T) dev.off()

setwd(paste0(DIR,"/Data"))

}


##################################################################################################################################################################################


# Appendix B. Alternative explanatory variables (UK 10-Year Bonds)


for (j in 1:2) {

if (j == 1) bonds_spread <- T
if (j == 2) bonds_spread <- F

setwd(paste0(DIR,"/Data")) 
bonds1 <- read.table("UK_bonds.txt", header = T, stringsAsFactors = F)
colnames(bonds1)[2:6] <- paste0(colnames(bonds1)[2:6], "_uk")
bonds2 <- read.table("DE_bonds.txt", header = T, stringsAsFactors = F)
colnames(bonds2)[2:6] <- paste0(colnames(bonds2)[2:6], "_de")
bonds3 <- read.table("FR_bonds.txt", header = T, stringsAsFactors = F)
colnames(bonds3)[2:6] <- paste0(colnames(bonds3)[2:6], "_fr")
bonds4 <- read.table("IT_bonds.txt", header = T, stringsAsFactors = F)
colnames(bonds4)[2:6] <- paste0(colnames(bonds4)[2:6], "_it")

bonds5 <- merge(bonds1[,1:2], bonds2[,1:2], all = T)
bonds5 <- merge(bonds5, bonds3[,1:2], all = T)
bonds5 <- merge(bonds5, bonds4[,1:2], all = T)
bonds5$Date <- as.Date(bonds5$Date)

bonds5$Price_mean <- (bonds5$Price_de + bonds5$Price_de + bonds5$Price_de)/3

for (cols in 2:ncol(bonds5)) {
  bonds5[,cols] <- na.approx(bonds5[,cols], n = length(bonds5[,cols]), na.rm = F, rule = 2)
}

bonds5$Price_uk_mean <- bonds5$Price_uk-bonds5$Price_mean

bonds <- bonds5
bonds <- subset(bonds, Date >= as.Date("2014-01-01")) 
bonds <- subset(bonds, Date < as.Date("2022-01-01")) 

for (cols in 2:ncol(bonds)) {
  bonds[,cols] <- na.approx(bonds[,cols], n = length(bonds[,cols]), na.rm = F, rule = 2)
}
bonds$Date <- as.Date(bonds$Date)

setwd(paste0(DIR,"/Plots"))

if (Print == T) if (bonds_spread == F) jpeg("rhs_3_Bond_yields.jpeg", units = "cm", width = 17, height = 16, res = 400)
if (Print == T) if (bonds_spread == T) jpeg("rhs_3_Bond_spreads.jpeg", units = "cm", width = 17, height = 16, res = 400)

if (bonds_spread == F) {
  plot(as.Date(bonds$Date), bonds$Price_uk, type = "l", col = "darkblue", 
     xlab = "Timeline", ylab = "Yield (%)", main = "UK 10-Year Government Bonds") # blue
}
if (bonds_spread == T) {
  plot(as.Date(bonds$Date), bonds$Price_uk_mean, type = "l", col = "darkblue", 
       xlab = "Timeline", ylab = "Yield (%) - Difference UK-EU Basket", main = "UK 10-Year Government Bonds") # blue
}

abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win elections => referendum coming
abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit

pos <- 3.03 
if (bonds_spread == T) pos <- 1.63
text(y = pos, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
text(y = pos, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
text(y = pos, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 

if (Print == T) dev.off()

setwd(paste0(DIR,"/Data"))

}


##################################################################################################################################################################################


### PLOTS II: LHS - Party Demands & Statements (Trends & Densities)


# Trends - Nationalist party communication

for (j in c(1,3)) {
  
  setwd(paste0(DIR,"/Data"))
  z <- read.table("LHS_plots_data.txt", header = T, stringsAsFactors = F)
  
  z$Date <- as.Date(z$Date)
  z2 <- z
  all_in_one <- T
  combined_dims <- F
  Density <- F
  
  setwd(paste0(DIR,"/Plots"))
  
  if (j == 1) {crit <- F; aggregate <- T}
  if (j == 2) {crit <- T; aggregate <- T}
  if (j == 3) {crit <- F; aggregate <- F}
  if (j == 4) {crit <- T; aggregate <- F}
 
  CNT <- c("AUT","DNK","FRA","DEU","ITA","NLD","FIN","POL","IRL","ESP")

  if (all_in_one == T) if (crit == F) if (aggregate == T) if (Print == T) jpeg("simple_aggr.jpeg", units = "cm", width = 17, height = 16, res = 400)
  if (all_in_one == T) if (crit == T) if (aggregate == T) if (Print == T) jpeg("simple_crit.jpeg", units = "cm", width = 17, height = 16, res = 400)
  if (all_in_one == T) if (crit == F) if (aggregate == F) if (Print == T) jpeg("all_in_one_aggr.jpeg", units = "cm", width = 25, height = 42.5, res = 400)
  if (all_in_one == T) if (crit == T) if (aggregate == F) if (Print == T) jpeg("all_in_one_crit.jpeg", units = "cm", width = 25, height = 42.5, res = 400)
  
  MAI <- c(.4,.6,.4,.1)
  if (aggregate == F) par(mfrow = c(5,2), mai = MAI)
  if (aggregate == T) {
    CNT <- NULL # for only aggregted plots
    options(warn = -1)
    par(mfrow = c(1,1), mai = MAI)
    par(default_params)
  }
  
  for (i in 1:length(CNT)) { 
    
    z <- z2
    if (i > 0) z <- subset(z, z$Country_Code == CNT[i])
    
    if (combined_dims == T) {
      z$Critique <- z$Critique/2
      z$Demands <- rowSums(cbind(z$Demands, z$Critique), na.rm=T)
      z$Demands[z$Demands == 0] <- NA
    }
    
    if (i == 0) CNT_NAME <- "Aggregate"
    if (i > 0) CNT_NAME <- unique(z$Country_Code)
    
    if (crit == F) {
      if (nrow(z[!is.na(z$Demands),]) > 5) {
        
        {
          SPAN <- .4
          if (i > 1) {
            if (CNT[i] %in% c("DEU","NLD")) SPAN <- .55
            if (CNT[i] %in% c("ESP","FRA","ITA","POL")) SPAN <- .7
            if (CNT[i] %in% c("DNK")) SPAN <- 1
            if (CNT[i] %in% c("FIN")) SPAN <- 10
          }
          l <- loess(z$Demands ~ as.numeric(z$Date), span = SPAN)
          newx <- seq(min(as.numeric(z$Date)), max(as.numeric(z$Date)), length.out=1000)
          pred <- predict(l, newx, se = T)
          fit <- pred$fit
          upr <- pred$fit+1.96*pred$se.fit
          lwr <- pred$fit-1.96*pred$se.fit
          preds <- cbind.data.frame(fit, lwr, upr)
          del <- which(is.na(preds$fit))
          
          NAME <- paste0(CNT_NAME, "_results.jpeg")
          if(all_in_one == F) if(Print == T) jpeg(NAME, units = "cm", width = 16, height = 16, res = 400)
          
          XLIM <- c(min(as.Date("2014-01-01")), as.Date("2022-01-01")+20)
          
          if(all_in_one == F) plot(jitter(as.numeric(z$Date), factor = 15), z$Demands, ylim = c(min(1,z$Demands,na.rm = T)-.2,4.35), xlim = XLIM, type = 'n', xlab = "", ylab = "Demands", xaxt = "n", yaxt = "n")
          if(all_in_one == T) plot(jitter(as.numeric(z$Date), factor = 15), z$Demands, ylim = c(min(1,z$Demands,na.rm = T)-.2,4.35), xlim = XLIM, type = 'n', xlab = "", ylab = "", xaxt = "n", yaxt = "n")
          title(main = paste0("Results: ", CNT_NAME))
          if(all_in_one == F) axis(1, at = as.Date(paste0(2014:2022,"-01-01")), labels = 2014:2022, cex.axis = 1)
          if(all_in_one == T) axis(1, at = as.Date(paste0(2014:2022,"-01-01")), labels = 2014:2022, cex.axis = 1)
          if (length(del) > 0) polygon(c(rev(newx[-del]), newx[-del]), c(rev(preds[-del,3]), preds[-del,2]), col = 'grey95', border = NA)
          if (length(del) ==0) polygon(c(rev(newx), newx), c(rev(preds[,3]), preds[,2]), col = 'grey95', border = NA)
          pal <- colorRampPalette(rev(c("red","orange","grey"))) #pal <- colorRampPalette(c("purple", "blue"))
          Col <- pal(nrow(z))[as.numeric(cut(z$Demands, breaks = nrow(z)))]
          set.seed(0)
          points(z$Date, jitter(z$Demands, factor = 1), pch = 19, cex = 1, col = scales::alpha(Col,.5))
          lines(newx, fit, col="black", lw=1.5) #lines(newx, upr, col="black", lw=1, lty = 2); lines(newx, lwr, col="black", lw=1, lty = 2)
          box()
          if(aggregate == T) axis(side = 2, at = c(1:4), labels = c("Status Quo", "Reform", "Reform or Leave", "Leave"))
          if(aggregate == F) axis(side = 2, at = c(1:4), labels = c(1:4), cex.axis = 1.3)
          
          abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
          abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win
          abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
          abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
          abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
          abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
          abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
          abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit

          text(y = 4.34, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
          text(y = 4.34, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
          text(y = 4.34, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 
          
          fit_integration <- fit
        }
        
        if(all_in_one == F) if(Print == T) dev.off()
        
      }
    }
    
    if (crit == T) {
      if (aggregate == T) z <- z2
      if (combined_dims == F) {
        {
          if (i != 7) {
            l <- loess(round(z$Critique) ~ as.numeric(z$Date), span = .5)
            newx <- seq(min(as.numeric(z$Date)), max(as.numeric(z$Date)), length.out=1000)
            pred <- predict(l, newx, se = T)
            fit <- pred$fit
            upr <- pred$fit+1.96*pred$se.fit
            lwr <- pred$fit-1.96*pred$se.fit
            preds <- cbind.data.frame(fit, lwr, upr)
            del <- which(is.na(preds$fit))
          }
          
          if(all_in_one == F) NAME <- paste0(CNT_NAME, "_criticism.jpeg")
          if(all_in_one == F) if (Print == T) jpeg(NAME, units = "cm", width = 17, height = 16, res = 400)
          
          plot(jitter(as.numeric(z$Date), factor = 15), round(z$Critique), ylim = c(.9,2.25),
               xlim = XLIM, type = 'n', xlab = "",
               ylab = "EU Criticism", main = paste0("Results: ", CNT_NAME), xaxt = "n", yaxt = 'n')
          axis(1, at = as.Date(paste0(2014:2022,"-01-01")), labels = 2014:2022, cex.axis = 1)
          axis(side = 2, at = c(1:2), labels = c("Criticism (policy)", "Criticism (institutional)"))
          if (i != 7) {
            polygon(c(rev(newx[-del]), newx[-del]), c(rev(preds[-del,3]), preds[-del,2]), col = 'grey95', border = NA)
            pal <- colorRampPalette(rev(c("darkblue","steelblue2"))) #pal <- colorRampPalette(c("purple", "blue"))
            Col <- pal(nrow(z))[as.numeric(cut(round(z$Critique), breaks = nrow(z)))]
            set.seed(0)
            lines(newx, fit, col="black", lw=1.5) #lines(newx, upr, col="black", lw=1, lty = 2); lines(newx, lwr, col="black", lw=1, lty = 2)
          }
          points(z$Date, jitter(round(z$Critique), factor = .5), pch = 19, cex = 1, col = Col)
          box()
          
          abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
          abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win
          abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
          abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
          abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
          abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
          abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
          abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit
          
          text(y = 2.22, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
          text(y = 2.22, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
          text(y = 2.22, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 
          
          if(all_in_one == F) if (Print == T) dev.off()
        }
      }
    }
    
  }
  
  if (Print == T) dev.off()
  
}


##################################################################################################################################################################################


# Densities:

{

z <- z2
CNT_NAME <- "Aggregate"
pal <- colorRampPalette(rev(c("red","orange","grey"))) #pal <- colorRampPalette(c("purple", "blue"))
rescale_by_freq <- F
require(scales)  

par(mfrow = c(2,2), mai = c(.4, .8, .4, .1))

Density <- F
if (Density == F) {
  
  for(i in 1:4) {
    
    xn <- z
    xn <- xn[!is.na(xn$Demands),]
    Shares <- table(xn$Demands)/sum(table(xn$Demands))
    shares <- Shares/max(Shares)
    #i=1
    xn <- xn[xn$Demands == i,]
    xn <- xn[order(xn$Date),]
    xd <- density(as.numeric(xn$Date), adjust = 1, bw = "sj")
    
    if (rescale_by_freq == F) rescale <- 1/max(xd$y)
    if (rescale_by_freq == T) rescale <- shares[i]/max(xd$y)
    
    if (rescale_by_freq == T) NAME <- paste0("ALL4_", i, "_density.jpeg")
    if (rescale_by_freq == F) NAME <- paste0("ALL4_scaled_", i, "_density.jpeg")
    if (Print == T) jpeg(NAME, units = "cm", width = 17, height = 16, res = 400)
    
    X_upper <- 1.05
    plot(xd$x, rescale*xd$y, ylim = c(-.1, X_upper), xlim = c(as.Date("2014-01-01"), max((z$Date))), 
         type = 'n', xlab = "", ylab = "Frequency of Demands (Density)", main = "", #yaxt = "n", 
         xaxt = "n")
    Col <- pal(13)[(4*i)-3]
    polygon(xd$x, rescale*xd$y, lty = 0, col = scales::alpha("gray95", .75))
    
    ydata <- (xd$y*rescale)
    ydata[ydata < .01] <- NA
    polygon(xd$x, ydata, lty = 1, lwd = 2, border = Col)
    abline(h = 0.005, lwd = 8, col = "white")
    box()
    
    title(main = paste0("Results: ", CNT_NAME))
    axis(1, at = as.Date(paste0(2014:2022,"-01-01")), labels = 2014:2022, cex.axis = 1)
    
    if (rescale_by_freq == F) {
      pc_of_total <- paste0(round(100*Shares[i], 1), " % of total\ndemands")
      text(y = .7, x = (as.Date("2014-07-01")+20), pc_of_total, col = "black", adj = 0, cex = 1)
    }
    
    abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
    abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win
    abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
    abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
    abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
    abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
    abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
    abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit
    
    text(y = X_upper, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
    text(y = X_upper, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
    text(y = X_upper, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 
    
    if (Print == T) dev.off() 
  }
}

}


# Relative Densities

{
  par(default_params)
  
  rd <- z2
  rd <- rd[!is.na(rd$Demands),]
  rd_lo <- rd[rd$Demands <= 2,]
  rd_hi <- rd[rd$Demands >= 3,]
  
  rd_lo <- rd_lo[order(rd_lo$Date),]
  D_lo <- density(as.numeric(rd_lo$Date), adjust = 1, bw = "sj")
  D_lo_d <- cbind.data.frame(as.Date(D_lo$x), D_lo$y)
  colnames(D_lo_d) <- c("Date", "Density_lo")
  rescale <- 1/max(D_lo_d$Density_lo)
  D_lo_d$Density_lo <- D_lo_d$Density_lo*rescale
  interp <- approx(x=D_lo_d$Date,y=D_lo_d$Density_lo, xout = seq(as.Date(min(D_lo_d$Date)), as.Date(max(D_lo_d$Date)), "days"))
  D_lo_d <- cbind.data.frame(as.Date(round(interp$x)), interp$y)
  colnames(D_lo_d) <- c("Date", "Density_lo")

  rd_hi <- rd_hi[order(rd_hi$Date),]
  D_hi <- density(as.numeric(rd_hi$Date), adjust = 1, bw = "sj")
  D_hi_d <- cbind.data.frame(as.Date(D_hi$x), D_hi$y)
  colnames(D_hi_d) <- c("Date", "Density_hi")
  rescale <- 1/max(D_hi_d$Density_hi)
  D_hi_d$Density_hi <- D_hi_d$Density_hi*rescale
  interp <- approx(x=D_hi_d$Date,y=D_hi_d$Density_hi, xout = seq(as.Date(min(D_hi_d$Date)), as.Date(max(D_hi_d$Date)), "days"))
  D_hi_d <- cbind.data.frame(as.Date(round(interp$x)), interp$y)
  colnames(D_hi_d) <- c("Date", "Density_hi")

  RD <- merge(D_lo_d, D_hi_d)
  
  setwd(paste0(DIR,"/Plots"))

  if (Print == T) jpeg("relative_densities.jpeg", units = "cm", width = 17, height = 16, res = 400)
  
  plot(RD$Date, c(RD$Density_hi-RD$Density_lo), type = "l", col = "orangered2", ylim = c(-.71, .3), lwd = 1.75,
       xlab = "", ylab = "Density aggresive - density moderate", main = "Relative Densities") # blue
  
  abline(v = as.Date("2013-01-23"), lty = 2, col = alpha("gray", .75)) # Referendum announced
  abline(v = as.Date("2015-05-07"), lty = 2, col = alpha("gray", .75)) # Tories win elections => referendum coming
  abline(v = as.Date("2016-06-23"), lty = 2, col = alpha("brown", .75)) # Brexit referendum
  abline(v = as.Date("2018-11-13"), lty = 2, col = alpha("gray", .75)) # First Withdrawal Agreement (Nov)
  abline(v = as.Date("2019-03-29"), lty = 2, col = alpha("gray", .75)) # Brexit Chaos
  abline(v = as.Date("2019-10-17"), lty = 2, col = alpha("gray", .75)) # Revised Withdrawal Agreement agreed (Johnson)
  abline(v = as.Date("2020-02-01"), lty = 2, col = alpha("gray", .75)) # Brexit
  abline(v = as.Date("2020-12-31"), lty = 2, col = alpha("gray", .75)) # Brexit
  
  pos <- .3 
  text(y = pos, x = (as.Date("2013-01-23")-20), "Ref.\nannounced", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2015-05-07")-20), "Tories win election,\nEU concessions", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2016-06-23")-20), "Brexit\nref.", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2018-11-13")-20), "Withdrawal agreement\n", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2019-03-29")-20), "\n1st Brexit Date", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2019-10-17")-20), "Rev.\nagr.", col = "black", adj = 1, cex = label_size) 
  text(y = pos, x = (as.Date("2020-02-01")+20), "Brexit\ndate", col = "black", adj = 0, cex = label_size) 
  text(y = pos, x = (as.Date("2020-12-31")+20), "End of trans.\nperiod", col = "black", adj = 0, cex = label_size) 
  
  if (Print == T) dev.off()
  
}



# R and package Versions:

# sessionInfo()
# R version 3.6.1 (2019-07-05)
# Platform: x86_64-w64-mingw32/x64 (64-bit)
# Running under: Windows 10 x64 (build 19044)
# 
# Matrix products: default
# 
# locale:
# [1] LC_COLLATE=English_Switzerland.1252  LC_CTYPE=English_Switzerland.1252   
# [3] LC_MONETARY=English_Switzerland.1252 LC_NUMERIC=C                        
# [5] LC_TIME=English_Switzerland.1252    
# 
# attached base packages:
# [1] stats     graphics  grDevices utils     datasets  methods   base     
# 
# other attached packages:
# [1] zoo_1.8-6     ggplot2_3.3.2
# 
# loaded via a namespace (and not attached):
# [1] Rcpp_1.0.2       lattice_0.20-45  withr_2.1.2      dplyr_1.0.0      crayon_1.3.4     R6_2.4.0        
# [7] grid_3.6.1       lifecycle_0.2.0  gtable_0.3.0     magrittr_1.5     scales_1.0.0     pillar_1.4.4    
# [13] rlang_0.4.6      rstudioapi_0.11  generics_0.1.0   vctrs_0.3.1      ellipsis_0.3.0   tools_3.6.1     
# [19] glue_1.4.1       purrr_0.3.4      munsell_0.5.0    compiler_3.6.1   pkgconfig_2.0.3  colorspace_1.4-1
# [25] tidyselect_1.1.0 tibble_3.0.1 
